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The twelve islands of the Choofootse group have to 
be visited each day for mail and supplies by an air crew. 
The islands are normally toured in the order ABCD...LM, 
since that ordering provides the shortest trip. 


During a labor dispute, the crew finds that the 
ordering of their daily trip is not specified, so they 
decide to visit the twelve islands in the order that will 
make the longest possible trip. What is that order? 


PROBLEM 68 


From the point indicated by the star, the twelve 
islands have coordinates as follows: 


48 = 5) -50, -40 
38 at -46, - 2 
Pf -54 -49, 24 


1 -43 -24, 39 


The problem is somewhat similar to Problem D3, The 
Pushbutton Radio Problem, in Problems for Computer Solution 
(Gruenberger and Jaffray), except that the distances 
involved are not equal, and the permutations are circular. 
Nevertheless, the analysis of Problem D3 indicates that 
the maximum trip between the twelve islands is of the 
order of 77 units, where one unit is the average distance 
between two successive islands. 


See also the problem of O'Gara, The Mathematical Mailman, in 
Martin Gardner's Sixth Book of Mathematical Games, page 235. Mr. 
Gardner refers also to Problem No. 64 in Hugo Steinhaus! One Hundred 
Problems in Elementary Mathematics. Ta 
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Book Review... 


INTEGER W(15).Z,CW 
DIMENSION 023) 
COMMON O 
WEITE(6-91) 
91 FORMAT (3UGH1IFORTRAN? A DISCOVERY APPROACH? 
al 17H BY ROBERT TEAGUES 
1 26H CANFIELD PRESS 1974// 
1 22H REVIEWED BY WeCeMCGcc//) 
CALL MPC 495e0,L,6) 
Ow=6 
N=C 
16 PEAD (5,92) (W(I),I=1,15) 
92 FORMAT (1515) 
OO 80 I=1,15 
Z=WC(I) 
IF €Z) 20,100,390 
20 Z=-Z 
N=1 
3 NCH=M00(7Z,10) 
NOW=Z/1¢ 
IF C(NOWLEQ.9) GO TO 4d 
CALL MP(241#MOD (CNOW-1,12) ptt (NOWH1LIS12_ NCH, OK+L) 
OW=OW+NCH 
GC TO 5) 
49 CALL HP(3,;NCHsU,L,OWH1) 
OW=OWF1 
50 IF«(NeEQ.1)? GO TD 60 
WRITE (6493) (O(J),J=1,0W) 
93 FCRMAT (1H »23A1) 
OwW=C¢ 
63 N=t 
80 CONTINUc 
GG TO 19 
itd STOF 
ENO 
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SUBROUTINE MP(SW,IyJoNyK) 
INTEGER SWeddeRqC 
DIMENSION TEXT(12.12) ,OUT(29) »PUNC (4) 
COMMON OUT 
DATA PUNC /4iHegiHsy91H',1H-/ 
BICL, Jd (1+ (dei) *229712¢1 
GO TO (1052u93t), SW 
19 REAG (5,11) TEXT 
41 FORMAT (12A1) 
RETURN 
29 R=I 
C=J 
00 25 L=1,N 
M=K+Lo4 
OUT(MID=TEXT (RC) (Data for this 
C=BU(R.C) Fortran program 
R=MOD(P,12) 41 is on the next 
25 CONTINUE page.) 
RETURN 
39 OUT(K) =PUNC(I) 
RETURN 
END 
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VENOTSHICONS 


ARESEDBOOKEA 
RKPIT 
Retencceces Mr. McGee's Fortran program outputs a review 
TOFASSPOSTSE 
SOCOMPROGFAM of Robert Teague's Fortran text. As one 
MINGODCETHER 
ATIVEXWHINAU might suppose, the review is quite favorable. 
TUDEMONMANGI 
CICALANGUAGE 
FORTRANOPLIU 
SISLLYTYWO3Y 
813 -762 1122 Si2 -812 52 194 520 492 $494) «6364 7634 -952 -742 2 
“1322 745 813 1217 872 -E€68 743° 1138 522 313 -e271 217221261 -7ed-il72 
1123 -12 “72 -453 4.1261 -592 -732 -4 ©5013 -843 led2 e422 -342 325 
S2-132e 172 (bye OBL SY tela e3 e121 -233 -382 743-343 2-572 -3U3"1521 
3 133 -553 172 1263-1292 -272 1113 -572 -+38c 473 2-142 -542 -153 152 
L2ge107e -13 al 813 1138 -892 -553 285 52 2-32 13 -294-1ic2 2 
“1201 13 492 §13 -463-1292 -851 1441 5ued 851 -633-1232 -8ie2 ACL obi 
Anas 2 9) Ges 2a BUA. eyelieolelaal PeSicke Sao La 2:52 <5 eee ES) ol aL 
=43 =451—1171 -S32 343 502 -533 1013 -668 745 -273 =-512-1362-1351 2 
eee SGyes » Biya 152 -?t2 772 1217 -005 -€38 ee Obie San elie 133-2243 1461 
“1)84-1u42-10 72-1061 855 -892 -3832-1492 593 1263 -663-1352-1ts2 133: 742 
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3-452 -233 1272 -892 -792 -172 133i 813 Ss) fal ose 3 Se 2-ieel 
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Speaking of Languages 


The column this month is short due to space 
restrictions. It would, however, be a shame to 
let the above book review go by without comment. 
The fact that the text and this column have the 
same author is not what is noteworthy, but rather 
the unusual form of the review. 

The Fortran program was written for IBM 
equipment and uses some non-standard techniques to 
produce the end result. This use of Fortran 
would appear to be unique, which might justify 
the use of these features. The program consists 
mainly of input/output statements, which read 
data cards and produce a printed listing that 
comprises the review. Take a close look at how 
the alphabetic material is produced; there is 
something to be learned from this exercise. @ 

All in all, this is not only a novel, but 
highly appropriate and educational means of 
producing a book review. 


ROBERT TEAGUE 


@ _ ART -: CONAPTIN Gear 


square 


MORE THAN YOU EVER WANTED To KNOW ABOUT root 


Square root, the most elementary of functions, can 
be calculated by many algorithms. (In what follows, the 
square root of 153 will be used as an example; its value is 


12 .3693168768529816494642295679222...) 


1. ‘The binomial theorem method. This is the 
algorithm that has been taught in our schools for hundreds 
of years, causing millions of students to hate arithmetic, 
and fear square root. The process is shown in Figure A. 


The number whose root is sought is divided into 
periods of two digits, both ways from the decimal point, 
For the leading period, a trial square root is taken 
(which requires knowing the squares up to 100). The 
trial root is squared and subtracted from the first period. 
The root up to this point is now doubled and used as a 
trial divisor. From here on, the process somewhat 
resembles long division. In the example, 20 is to be 
divided into 53, but with the restriction that the quotient 
is to be first added to the trial divisor. Thus, 53 
divided by 20 is 2, but the trial divisor changes to 
20 + 2, as shown, New periods are brought down as needed, 
and the process repeats, getting progressively more awkward. 


The geometric analog is shown, at the stage where 
the root has been located up to the decimal point, in 
Figure B. The large square has an area of 153 units, 
an@ the square of 12 (crosshatched) can be removed. What 
is left consists of two parts. We are expressing 153 as 
& the square of something: 


153 = (12 + x)® = Wh + 2126x + x? ® 
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and the x* term is disregarded as momentarily unimportant. 
. The doubling in the algorithm is the 2 of the binomial 
expansion. We seek the width of the strips that represent 
2°12-x. When we subtract 144 from 153, the difference of 
9 is made into 900 by bringing down another period, and we 
"divide" by 24 to determine the width. As with many other 
methods, at all stages we really deal with the sequence of 
the digits involved; that is, we turn the problem into 
a problem in integers at all times. 


The process is awkward and inefficient and requires 
memorizing an abnormal number of arbitrary rules. Since 
most seventh grade teachers have no idea why it works (but 
who recall that they too hate it), this fear and loathing 
is readily communicated to the students. And this is all 
quite mysterious, inasmuch as the algorithm given below 
has been around for hundreds of years and is superior in 
every way. 


2. Newton's Method. This algorithm, which is the 


Newton-Raphson scheme for finding the root of an equation 
(applied to x© - N = 0) is usually given as an iterative 
formula: 


N 
Xsay = he Nees + x, ) 


The method predates Newton by many years. It can be 
explained in these terms: 


a 


» fw 


Using division, the number N is divided by any trial 


i 
divisor, p, yielding a quotient q. i aa oe -Wae (rnd ol =) oon Z A 
the definition of square root, p is the desired root. nN 
In general, p will not equal q, and therefore one of p fe 
or q is less than the square root and the other is greater; 
that is, the root has been located between p and q. Any 
value between p and q will be a better approximation to the 
root and, lacking any other information, we take (ptq)/2 
as the next approximation (what Newton showed, in fact, 
was that (p+q)/2 is the best choice that can be made simply.) 

The process now repeats, using (p+q)/2 as the new trial 
divisor. The arithmetic for N = 153 is shown in Figure OC, 
using 12 as a starting value. 

1275 

wey] 153 1200 + 1275 

144. a LT ae S238 
90 e 
B4 
60 
123586 
1238 | 153000000 
1238 C 
2920 123800 + 123586 
2476 ne. Men, = 123693 
LALO 
3714 
7260 Of the two possible sequences for the 
6190 Square root of the given sequence, we 
10700 have established the one we want; hence, } 
9904 the decimal points are of little 
consequence in the calculation. 
7960 
7428 
5320 


The Newton algorithm has the following advantages: 


a. The process is fast. Stated in crude terms, the 
method allows us to double the number of correct digits at 
each stage. For each division, the quotient should be 
calculated to twice the number of digits that repeat. For 
example, in Figure C, the digits 12 repeat in the first 
division, and hence the quotient is carried to 4 digits. 
The result is not quite good to 4 digits, as seen in the 
second stage, where the digits 123 repeat. The second 
quotient is carried to 6 digits and, after interpolation, 
the result should be good to about 6 digits, which it is. 


self-correcting. If an error is made, then the result at 


b. Since the procedure is iterative, it 1s also » 
that stage is simply not as good as it might be. 


| 
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ec. The procedure stems directly from the definition . 
of square root, and hence has few rules to remember. 


dad. The scheme capitalizes on previous information. € 
For example, if it is known that the square root of 2 
begins with 1414, then that value can be used as the 
starting value, and one division (and interpolation) will 
yield 8 digits of the root. Thus, using a mechanical 
desk calculator and a one-page table of roots good to 4 
significant digits, any square root can be obtained to 
the machine's limit (10 digits) with two divisions. And 
since the algorithm is based on division, it utilizes the 
part of the desk calculator for which the most money was 
spent. 


e. The method is readily programmed for any computer, 
in any language. To avoid the problems of judgement 
invoived in selecting a suitable starting value and in 
testing for convergence of the process, it is expedient to 
use either one or N as the starting value (which really 
makes the initial value (1+N)/2 and taking a fixed number 
of iterations, with each division carried out to the 
precision required of the result. Thus, in Fortran, 
using 8-digit floating arithmetic, six iterations will 
handle the worst case, using the starting value above. 


3.. Other approximation methods. Any of the weli- 
known schemes for locating a root of an equation can be 
applied to the equation x° - N= 0. Horner's method is 
out of style, but the false position method is still used. 
See Figures D and E, We attempt to locate the point at 
which the curve y = xe 153 crosses the x-axis. The 
region between x = 12 and x = 13 has been greatly enlarged. 
On the assumption that the straight (dotted) line is an 
approximation to the curve, the similar triangles involved 
imply that 


Z = Nie) 
Xo - Ki 7 Yot yl 
where all the values are taken as positive. We have then 


Xo - %, = 1, Yo= ia o i155} = 16.5 Vawer Los 12° = 9. 
From the relation above, Z = 16/25, or .64, and our next 
approximation is then 12.36. It is not at all clear to 
what extent we have located the root, but we can calculate 
a few values of the function, and determine that the 
function changes sign between 12.36 and 12.37. The 
process can now repeat, using those values as x, and Xo. 


The false position method seems best suited to hand 
calculation. 


4, Logarithms. If we divide the logarithm of 153 
by 2: 


log 153 = 2.18469143081759880313022... 
log Y153 = 1.09234571540879940156011... 


the antilog of the result, found by interpolating in a 
7-place table of logarithms (Vega's table) is 12.3693. 


5. Hastings Approximations. Both the logarithm 
and the antilogarithm can be obtained from the ingenious 
rational functions found in Cecil Hastings! book, 
Approximations for Digital Computers (Princeton University 
Press, 1955). The logarithm given in (4) above was 
determined by direct (series) calculation of the logarithm. 
Hastings gives (his Sheet 20) the following approximation: 


for O<x<1l 


oe 2 Slot 7° 
10 GPE aa et tao X okt gx treet alx ) 
where a. = 1.15129277603 
1 a5 = .01742111988 
a, = .66273088429 ash ete 
a3 = .25439357464 poe oe 
ay = 07295173666 a7 = Ce 7 
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We have here x = .09234571541. By substituting into 
Hastings! formula, the square root of 153 is found to be 
12.3693168. 


6. Table lookup. This is the first algorithm that 


most people meet; namely, to find a table of square roots 
and look up the desired answer. If a table of roots is 
stored in a computer, the table can be searched (linearly, 
or by binary searching), or the address of the required 
function value can be itself calculated. For example, 

if a table of 1000 entries (000 through 999 for the 
argument) is stored in a block addressed at T, then the 
functional value for 153 is at address T4153. 


7. Interpolation. Given a table of square roots 


ana differences, as follows: 


N SQRT Successive differences 
153 123693169 (4) 

+ 
154 124096736 403567 


(=) 
155 124498996 402260 1307 


156 124899960 4.00964 1296 alah 

157 =: 125299641 399681 1283 13 

158 125698051 398410 1271 12 
we can find the square root of other values within the 
game range by utilizing well-known schemes for interpolation. 
To find the square root of 153.5, for example, Newton's 
forward interpolation formula specifies that 

X - Xo A5}5}515) co allsys} 


h a eae ete aa 1/2 


3 
Ryotes: 


men) w= 1153.5 2 Yo t+ udyo + ae Avy = ee 


(0) 
where the a's indicate successive differences. We have here 
f(u) = 123693169 + (1/2)(403567) + (-1/8)(-1307) +--- 

= 123895114 


which is correct to 9 significant digits. 


8. Bracketing methods. in a computer, the range 
of numbers to be considered is always restricted. For 
example, if we are working in fixed point on a machine with 
a 24-bit word, then the possible range on N is from zero 
to 8388607, which means that any square root must lie 
between zero and 2896, A perfectly straightforward 
algorithm, then, would be a loop to test the square of 
every value between zero and 2896 to determine which one 
most closely yields N. This would be terribly inefficient. 
The procedure shown in Figure F will speed up the process 
by a crude form of interval-halving. The resulting 
algorithm is discovered independently by countless sixth 
graders who want square roots but have not yet been exposed 
to any algorithm, even the miserable one give in (1) above. 


As shown on the flowchart, the range on N has been 
further restricted to 20 bits, and the calculation for 
N = 15300 produces these numbers: 


LO TRY HI 
0 0 999 
) 499 999 
) 499 499 
0 249 ag 
0 124 12k 
0 62 = 124 
62 62 12k 
93 o3 ee 
108 108 = «124 
116 1160s 124 
120 120 124 
123 123° 9 tel: 


and when the process converges, the root has been located 
within the limits of the specified precision. The algorithm 
is quite inefficient, but has two advantages: 


a. it is quite easy to program. 


b. It uses no division. The "averaging" step 
can be done by multiplying by .5, or, in a binary computer, 
by a right shift of one bit position. 


If we deal with numbers in scientific notation in the 
computer, then the admissable range on N is much greater, 
and the preceding scheme would be hopelessly inefficient. 

A more efficient bracketing scheme is shown in Figure G. A 
search is conducted for the root, starting at zero and 
proceeding by large steps (D= 10). The test determines 
when the root has been passed. When this occurs, the 
procedure calls for backing off one full step, and the 
search increment is then reduced by a factor of 10. A 
counter, K, tallies the number of times this has taken place; 
each such occurence yields one digit of precision. 
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Housekeeping: Set L@ = 000 
TRY = 000 
Bie & Ses ; 
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TRY —> HI 


Set TRY = 


average of 
L@ and HI 


Continue 
(X is the root) 


The algorithm of Figure G is also easy to program, 
and also uses no divisions. It is shown on the flowchart 
in general form. In practice, the value of the search 
interval and the initial search point should be functions 
of N, to cut down the machine time. 


9. Odd number subtraction. From a table of 
differences of consecutive squares, or from the identity: 


MES GEIS era a ane. (2s Werle 


it can be seen that the sum of consecutive odd numbers 
equals the square of the consecutive integers. This fact 
can be used in reverse. By subtracting the odd numbers 
from a number whose square root is desired, and counting 


the number of subtractions, the root can be obtained. For 


example, for N = 153, successive subtractions of odd 
integers produce this pattern: 


start 153 count the subtractions: 
-1 152 al 
-3 149 2 
“5 144 3 
-7 Si 4 
-9 128 5 
-11 117 6 
-13 104 1 
-15 89 8 
-17 72 9 
-19 DS 10 
-21 32 11 
-23 9 12 
=25 underflow 


indicating that, to the nearest integer, the square root 
aig dle, To extend the process, one begins with 15300 or 


“ 1530000. The process works from left to right, similarly 


to long division. 


PC20-13 


Before the days of automatic division on desk 
calculators, the odd number subtraction method was widely 
used for square root, particularly on kKey-driven machines 


(e.g., the Comptometer). Around 1952, the Friden company 


tried to mechanize the process on the first square root 


desk calculator, but the result was a failure. 


The user 


of that machine would enter a number into the lower dial, 
using the dividend entry key. When the square root key 
was depressed, the machine would start subtracting successive 
odd numbers from the lower dial. When that dial rolled 
negative, the machine would reverse and add one cycle 
(just as in division), shift the carriage, and begin 
subtracting odd numbers again. The trouble with the 
procedure was that two banks of the keyboard mechanism 

had to be controlled, as in the example we have seen. The 
mechanical details of this procedure were complex, and the 


machine broke down rapidly and frequently. 


» 


, The solution to the problem, which appeared around 
1954 in a new Friden square rooter, was ingenious. When 
the square root key on the new machine was depressed, the 
machine first multiplied N by 5, and then proceeded to 
subtract odd multiples of 5. The process is mathematically e 
the same, but mechanically simpler, since the series of 
numbers 05, 15, 25, 35, 45,...,95 can be generated by 
locking a 5 in one bank of the keyboard and then generating 
consecutive integers in the bank to its left. The later 
machine worked fairly well, and for years was the only 
desk machine on the market with square root capability. 
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10. The quadratic method. The roots of the 
quadratic equation 


x° =x +M = 0 (A) 
can be found by the quadratic formula. If the roots can 


also be found by some other means, then, since the formula 
involves a square root, there results a very roundabout 


way of taking square roots. The equation (A) can be 
transformed into an interative scheme: 
cS x* + M 
aS (B) 
X54 = xT + M 
To find the square root of a number N, form the quantity 
l-WN 
iQ =—= 
4 


with N chosen in the range from zero to one (which is not 
a serious restriction, since we seek the sequence of digits 
for the root and can adjust the decimal point later). 


Thus, to find the square root of 153, take N = .0153. 

Then M= .246175. That value for M is used in equation (B), 
with some suitable starting value for x, and the process is 
iterated to convergence. When it converges, the required 
root will be (1 - 2x). The convergence rate is abysmally 
slow. The process seems useless as a practical method 

for finding square roots; it is included simply because it 

is a method, 


ll. Graphical methods. Figure H shows how any 
square root can be built up geometrically. The diagonal 
of a unit square is the square root of 2. If that length 
is brought down (are AB, centered at 0) and a perpendicular 
is erected, then the new diagonal is the square root of 3, 
and so on. Again, this method is somewhat impractical, 
but is another algorithm for square rooting. 
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12. The binomial theorem again. For numbers that 
differ from a perfect square by a square, the binomial 


theorem gives another method for square root. For example, 
1/2 = 
(ays 9)" = (anny? 4 (aye) (aun) /?(9) 
+ (1/2)(-1/2) J37/3 
emer co we 
(1/2) (-1/2)(-3/2) : 
mere ieamens Gl eS Bec 


2:3 


12 + .375 = .0058593 + .0001831 - .O0000071 ++++ 
12 .3693167 


13. Some of these schemes can be combined. For 


example, table lookup can be used to furnish a good 
starting value for any of the iterative methods. For 
many years, the makers of mechanical desk calculators 
furnished one-page tables of square roots, with directions 
for extending the precision by the Newton method. 


Every subroutine library and every programming 
package has a square root method built in. Which method 
is used is partly a matter of taste on the part of the 
author of the package. In a floating point system, it 
is most probably a Hastings approximation, perhaps followed 
by one stage of the Newton method. 


Users of computers are indebted to many brilliant 
men who have brought the art of computing to its present 
level. In particular, though, we all use the work of 
three men constantly: John von Neumann, who started it 
all; Jay Forrester, who invented the magnetic core; and 
Cecil Hastings, Jr., whose approximations are at work in 
every function subroutine package. a 


Printed by Argold Press, Inc. 
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«3010299956639811952137 388947 2449 3026768189881462109 
«9957 3227 355399099 34 3522357 614254077 5676601622989028 
.472135954999579 392818 347 3374625524 70881236719223051 
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e 485165195 .4097902779691068 3054 1540558684 6389889448472 


5435361080031597799614270974016597985065275 


7 8769956796 . 0826994747 5225559 3703897066064114447195437 


i 


tan 20 il 


———s 
——— 


24.34209842605184123904. 35479909902 349851867 


. 52083793107 2953857821 3154.04604906560607 307619264046 
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